function fout = DetectObjects(nstartframe, nendframe, nfindobjs, nfindmaxobjs,...
                              fvideo, fout, fbackground, fsubtract, handles, nDrawFrameSelect, nDrawEveryFrames)

%objs: contains detected objects for each frame. Is arranged such that: objs(frame number;
%object number; x1 x2 y1 y2), each column represents an x or y coordinate and is organized 
%as: x1 x2 y1 y2. Where x1, x2 are the vertical values (x1 = top, x2 = bottom) and y1 
%and y2 are the horizontal values (y1 = left, y2 = right ) of the square that depicts 
%the location of the object - this vert./horz.<->x/y association was selected
%due to the way matlab associates images and arrays. 

%areas: contains areas within the object (if applicable) - e.g. a beak on a bird or
%an area of the animal that is unique to male or female. This 3d-array is organized such 
%that: areas(frame number; object number; x1, x2, y1, y2, ntype, score). ntype
%represents the type of area: e.g. beak = 0, male area = 1. score provides
%a measure of the detection algorithm's likelihood of having detected the
%area successfully. 


load(fsubtract);
load(fbackground);

v = mmreader(fvideo);
cc = [cr1, cr2, cg1, cg2, cb1, cb2, ch1, ch2, cs1, cs2, cv1, cv2];

NumberOfFrames = get(v,'NumberOfFrames');

%if bupdatebackground equals 1 then the background is updated every nframesnewbg 
bupdatebackground = 0; 

%minimum number of pixels required for an object to be detected (otherwise
%considered noise)
nminpixels = 0; 


img = read(v,1);
imgsz = size(img);
imgx = imgsz(1);
imgy = imgsz(2);
imghalfx = round(imgx/2.0);
imghalfy = round(imgy/2.0);

xl = imgx;
yl = imgy;
r = 255;
g = 255;
b = 255;
linethickness = 3;

nstart = nstartframe; 
nend = nendframe;

if nendframe == -1
    nend = NumberOfFrames; 
end

nframes = nend-nstart + 1;


vx1s = zeros(nframes, nfindmaxobjs);
vx2s = zeros(nframes, nfindmaxobjs);
vy1s = zeros(nframes, nfindmaxobjs);
vy2s = zeros(nframes, nfindmaxobjs);

vex0s = zeros(nframes, nfindmaxobjs);
vey0s = zeros(nframes, nfindmaxobjs);
vephis = zeros(nframes, nfindmaxobjs);
veas = zeros(nframes, nfindmaxobjs);
vebs = zeros(nframes, nfindmaxobjs);
vets = zeros(nframes, nfindmaxobjs);

objs = zeros(nframes, 1);
vxheights = zeros(nframes, nfindmaxobjs);
vyheights = zeros(nframes, nfindmaxobjs);
vareas = zeros(nframes, nfindmaxobjs);


vcount = 0;

stepcount = 0;

totalframes = nend-nstart+1;


%Any areas that should be exluded from detection (x1a to x2a, y1a to y2a, 
%x1b to x2b, and y1b to y2b)
excludex1a = 0;
excludex2a = 0;

excludex1b = 0;
excludex2b = 0;

excludey1a = 0;
excludey2a = 0;

excludey1b = 0;
excludey2b = 0;



dimg0 = double(backgroundimg);
dhimg0 = double(rgb2hsv(backgroundimg));


%-If bupdatebackground is equal to 1 (set at the beginning) a new 
%background image is calculated every nframesnewbg, and imbuf and 
%countbuf are arrays used to calculate the new background image.
%-The background image is caluclated by adding frames every
%nframesupdatebuf. Objects are detected in each frame and only the
%background is detected. For each object add nbuf pixels to help ensure
%no foreground pixels are included. 
%-When a new image is caluclated we check that is likely to be correct by
%sutracting it from the old image, convert the subtract image to binary, 
%and count the number of white pixels. If the number of white pixels is
%less than nmaxbwpixels we keep the image (otherwise use the old image
%until the next update). 
imbuf = zeros(imgx,imgy,3);
countbuf = zeros(imgx,imgy);
nbuf = 10;
nframesnewbg = 300; 
nframesupdatebuf = 10;
nframecount = 0;
nmaxbwpixels = 10; 


for i = nstart:nend 
    
    if bupdatebackground == 1 && mod(i, nframesnewbg) == 0 && i > nstart
        
        ncount = find(countbuf(:) == 0, 1, 'first');        
        
        if isempty(ncount)
        
            imbuf(:,:,1) = imbuf(:,:,1)./countbuf;
            imbuf(:,:,2) = imbuf(:,:,2)./countbuf;
            imbuf(:,:,3) = imbuf(:,:,3)./countbuf;

            newbackgroundimg = uint8(imbuf);        
        
                
            [bwimg, idif] = ProcessImage(newbackgroundimg, dimg0, dhimg0, cc, nstremove, nstfill, nbwarea, dthresh); 
                
            nbwpixels = sum(sum(bwimg));
        
        
            if nbwpixels < nmaxbwpixels
                imgupdated = i
                dimg0 = double(newbackgroundimg);
                dhimg0 = double(rgb2hsv(newbackgroundimg));
            else
                imgnotupdated = i
            end
                
        end
        
        imbuf = zeros(imgx,imgy,3);
        countbuf = zeros(imgx,imgy);        
        nframecount = 0;
    end
        
        
    img = read(v,i);    %current frame's imgage

    %subtract background image (dimg0) from current frame's image
    [bwimg, idif] = ProcessImage(img, dimg0, dhimg0, cc, nstremove, nstfill, nbwarea, dthresh);  
        
    
    %remove any areas that have been specified for removal
    if excludex1a > 0 & excludex2a > excludex1a
        bwimg(excludex1a:excludex2a,:) = 0;
    end
    
    if excludex1b > 0 & excludex2b > excludex1b
        bwimg(excludex1b:excludex2b,:) = 0;
    end
    
    if excludey1a > 0 & excludey2a > excludey1a
        bwimg(:, excludey1a:excludey2a) = 0;
    end
    
    if excludey1b > 0 & excludey2b > excludey1b
        bwimg(:, excludey1b:excludey2b) = 0;
    end
    
    
    %locate objects in binary image
    CC = bwconncomp(bwimg);
    L = labelmatrix(CC);
    s  = regionprops(L, 'centroid', 'Orientation', 'MinorAxisLength', 'MajorAxisLength', 'FilledArea');
    
    cd = cat(1, s.Centroid);
    ori = cat(1,s.Orientation);
    mil = cat(1,s.MinorAxisLength);
    mal = cat(1,s.MajorAxisLength);           
    npixels = cat(1,s.FilledArea);               
    
    nobjs = max(max(L));

    img2 = img;    
    
    x1s = zeros(nobjs,1);
    x2s = zeros(nobjs,1);
    y1s = zeros(nobjs,1);
    y2s = zeros(nobjs,1);
    
    ex0s = zeros(nobjs,1);
    ey0s = zeros(nobjs,1);
    ephis = zeros(nobjs,1);
    eas = zeros(nobjs,1);
    ebs = zeros(nobjs,1);
    ets = zeros(nobjs,1);
    
    ncount = 0;
    
    areas = zeros(nobjs,1);    
    keepobjs = 1:nobjs;
    
    fx1s = zeros(nobjs,1);
    fx2s = zeros(nobjs,1);
    fy1s = zeros(nobjs,1);
    fy2s = zeros(nobjs,1);
    
    for j = 1:nobjs
        [row,col] = find(L == j);
        [x1, indx1] = min(row);
        [x2, indx2] = max(row);
        [y1, indy1] = min(col);
        [y2, indy2] = max(col);
        
        fx1s(j) = x1;
        fx2s(j) = x2;
        fy1s(j) = y1;
        fy2s(j) = y2;
        
        xheight = x2 - x1;
        yheight = y2 - y1;
        area = xheight*yheight;
        areas(j) = area;        
    end

    if nobjs > nfindmaxobjs
        keepobjs = zeros(nfindmaxobjs, 1);
        tempareas = areas;
        nk = 0;
        
        for k = 1:nfindmaxobjs
            [m, mindex] = max(areas);
            
            if mindex > 0
                areas(mindex) = 0;
                nk = nk + 1;
                keepobjs(nk) = mindex;
            end
        end
        
        areas = tempareas;
        
    end    
 

    for j = 1:nobjs             
        x1 = fx1s(j);
        x2 = fx2s(j);
        y1 = fy1s(j);
        y2 = fy2s(j);
        area = areas(j);      
        
        if npixels(j) > nminpixels & ~isempty(find(keepobjs == j, 1, 'first'))
            ncount = ncount + 1;
            
            x1s(ncount) = x1;
            x2s(ncount) = x2;
            y1s(ncount) = y1;
            y2s(ncount) = y2;                                            
            
            ex0s(ncount) = cd(j, 1);
            ey0s(ncount) = cd(j, 2);
            
            t = ori(j);
            phi = (abs(t-90)/180)*pi;
            
            ets(ncount) = t;
            ephis(ncount) = phi;
            
            eas(ncount) = mil(j)/2;
            ebs(ncount) = mal(j)/2;            
        end
    end    
  
    if ncount < nobjs        
        x1s(ncount+1:nobjs) = [];
        x2s(ncount+1:nobjs) = [];
        y1s(ncount+1:nobjs) = [];
        y2s(ncount+1:nobjs) = [];    
        
        ex0s(ncount+1:nobjs) = [];
        ey0s(ncount+1:nobjs) = [];
        ephis(ncount+1:nobjs) = [];
        eas(ncount+1:nobjs) = [];
        ebs(ncount+1:nobjs) = [];
        ets(ncount+1:nobjs) = [];
    
        nobjs = ncount;
    end
  
   
    if nobjs > nfindmaxobjs 
     nobjs = nfindmaxobjs;
    end
    
    vcount = vcount + 1;
    
    objs(vcount) = nobjs;
    
    for j = 1:nobjs         
        
        x1 = x1s(j);
        x2 = x2s(j);
        y1 = y1s(j);
        y2 = y2s(j);
        
        xheight = (x2-x1);
        yheight = (y2-y1);
        
        vx1s(vcount, j) = x1;
        vx2s(vcount, j) = x2;
        vy1s(vcount, j) = y1;
        vy2s(vcount, j) = y2;
        
        vex0s(vcount, j) = ex0s(j);
        vey0s(vcount, j) = ey0s(j);
        veas(vcount, j) = eas(j);
        vebs(vcount, j) = ebs(j);
        vephis(vcount, j) = ephis(j);
        vets(vcount, j) = ets(j);                
        
        vxheights(vcount,j) = xheight;
        vyheights(vcount,j) = yheight;
        vareas(vcount,j) = xheight*yheight;
                    
        
        if nDrawFrameSelect > 0                       
            if nDrawFrameSelect == 1
                img2 = DrawThickColorRGBBox(img2,  x1, y1, x2, y2, r, g, b, linethickness, xl, yl);
            else
                x0 = ex0s(j);
                y0 = ey0s(j);
                phi = ephis(j);
                a = eas(j);
                b = ebs(j);
                img2 = DrawImgEllipse(img2, x0, y0, a, b, phi);                            
            end
        end        

    end
           
    %If we are updating the background
    if bupdatebackground == 1 && nframecount > nframesupdatebuf && nobjs == nfindobjs
               
        bwtemp = zeros(imgx,imgy);                        

        for j = 1:nobjs                 
            x1 = x1s(j);
            x2 = x2s(j);
            y1 = y1s(j);
            y2 = y2s(j);
            
            x1 = x1-nbuf;
            x2 = x2-nbuf;
            y1 = y1-nbuf;
            y2 = y2-nbuf;
            
            if x1 < 1 
                x1 = 1;
            end
            
            if x2 > imgx 
                x2 = imgx;
            end
            
            if y1 < 1 
                y1 = 1;
            end
            
            if y2 > imgy 
                y2 = imgy;
            end
            
            bwtemp(x1:x2,y1:y2) = 1;
        end
         
        dimg = double(img);
        
        for j = 1:imgx
            for k = 1:imgy                
                if bwtemp(j,k) == 0                    
                    imbuf(j,k,1) = imbuf(j,k,1) + dimg(j,k,1);
                    imbuf(j,k,2) = imbuf(j,k,2) + dimg(j,k,2);
                    imbuf(j,k,3) = imbuf(j,k,3) + dimg(j,k,3);
                
                    countbuf(j,k) = countbuf(j,k) + 1;
                end
            end
        end         
        
        nframecount = -1;
         
    end

    nframecount = nframecount + 1;
    
     if nDrawFrameSelect > 0        
        if mod(stepcount, nDrawEveryFrames) == 0
            imshow(img2)
        end
     end
    
    stepcount = stepcount + 1;
    p = floor(100*(stepcount/totalframes));
    s = strcat('Processing Video...', num2str(p));
    s = strcat(s, '%');    
    set(handles.textoutput, 'String', s);        
    
end


nprocessstart = nstart;

save(fout, 'objs', 'vx1s', 'vx2s', 'vy1s', 'vy2s', 'vex0s', 'vey0s', 'vephis', 'veas', 'vebs', 'vets', ...
    'vxheights', 'vyheights', 'vareas', 'nprocessstart', 'imgx', 'imgy', 'imghalfx', 'imghalfy', ...
    'nfindobjs', 'nfindmaxobjs', 'nstartframe', 'nendframe');



